{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Using PyMC3 samplers on PyMC4 models"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pymc4 as pm\n",
    "import tensorflow as tf\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create simple pymc4 model\n",
    "@pm.model(auto_name=True)\n",
    "def t_test():\n",
    "    mu = pm.Normal(0, 1)\n",
    "\n",
    "model = t_test.configure()\n",
    "\n",
    "model._forward_context.vars\n",
    "func = model.make_log_prob_function()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create function to evaluate logp and dlogp over array of inputs\n",
    "@tf.function\n",
    "def logp_array(array):\n",
    "    #mu = array[0]\n",
    "    with tf.GradientTape() as tape:\n",
    "        tape.watch(array)\n",
    "        logp = func(array)\n",
    "    grad = tape.gradient(logp, array)\n",
    "    \n",
    "    return logp, grad"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# As the above function expects TF inputs and outputs, wrap it as PyMC3's samplers want numpy\n",
    "def logp_wrapper(array):\n",
    "    logp, grad = logp_array(tf.convert_to_tensor(array))\n",
    "    return logp.numpy(), grad.numpy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pymc4.hmc import HamiltonianMC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "size = 1\n",
    "n_samples = 500"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "tf.random.set_seed(123)\n",
    "np.random.seed(123)\n",
    "hmc = HamiltonianMC(logp_dlogp_func=logp_wrapper, size=size, adapt_step_size=False)\n",
    "curr = np.ones(size, dtype='float32') * .05\n",
    "posterior_samples = []\n",
    "stats = []"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "0.25\n",
      "10\n",
      "2.4338525715809687\n",
      "20\n",
      "0.8835782022550135\n",
      "30\n",
      "2.1327029448576513\n",
      "40\n",
      "0.6636198591746387\n",
      "50\n",
      "0.7703623177523709\n",
      "60\n",
      "1.8082085533140866\n",
      "70\n",
      "2.016802124717785\n",
      "80\n",
      "0.6217015104339201\n",
      "90\n",
      "1.207940992513313\n",
      "100\n",
      "1.3240099603905968\n",
      "110\n",
      "1.346048853571611\n",
      "120\n",
      "1.4427577775129503\n",
      "130\n",
      "0.8565491707334623\n",
      "140\n",
      "0.9624691955794613\n",
      "150\n",
      "2.4720718932582906\n",
      "160\n",
      "1.5805445755750056\n",
      "170\n",
      "2.1380879550318665\n",
      "180\n",
      "1.420551901927288\n",
      "190\n",
      "1.772666882095104\n",
      "200\n",
      "1.4909935183397898\n",
      "210\n",
      "0.9564664645510188\n",
      "220\n",
      "1.0094633430712847\n",
      "230\n",
      "0.9415898885107509\n",
      "240\n",
      "0.9789093380505235\n",
      "250\n",
      "1.006542783540078\n",
      "260\n",
      "0.9103849564026153\n",
      "270\n",
      "1.55622156748228\n",
      "280\n",
      "3.088911988511081\n",
      "290\n",
      "1.4383936543460296\n",
      "300\n",
      "1.1121804701838218\n",
      "310\n",
      "1.424086887611018\n",
      "320\n",
      "2.20499869286127\n",
      "330\n",
      "1.519518011454376\n",
      "340\n",
      "1.3419589165701242\n",
      "350\n",
      "1.9586890867688238\n",
      "360\n",
      "1.6645487485996129\n",
      "370\n",
      "0.9915868782338413\n",
      "380\n",
      "1.0335682127160022\n",
      "390\n",
      "2.580381765734188\n",
      "400\n",
      "0.8539483762542414\n",
      "410\n",
      "1.6204564049940595\n",
      "420\n",
      "1.7482972684183884\n",
      "430\n",
      "1.0932177164504921\n",
      "440\n",
      "1.1029667164138135\n",
      "450\n",
      "1.2918166766210588\n",
      "460\n",
      "1.2928383047613468\n",
      "470\n",
      "1.6941335213881863\n",
      "480\n",
      "2.3599326552982838\n",
      "490\n",
      "1.3764049670244904\n"
     ]
    }
   ],
   "source": [
    "# %%time  # NB: uncommenting cell magic %%time will prevent variable from escaping local cell scope\n",
    "\n",
    "for i in range(n_samples):\n",
    "    curr, stat = hmc.step(curr)\n",
    "    posterior_samples.append(curr)\n",
    "    stats.append(stat)\n",
    "    if i % 10 == 0:\n",
    "        print(i)\n",
    "        print(hmc.step_size)\n",
    "    \n",
    "trace = np.array(posterior_samples)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compare with `PyMC3`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "WARNING: Logging before flag parsing goes to stderr.\n",
      "W0414 17:16:25.698212 140348081661760 configdefaults.py:1458] install mkl with `conda install mkl-service`: No module named 'mkl'\n"
     ]
    }
   ],
   "source": [
    "import pymc3 as pm3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "with pm3.Model() as model3:\n",
    "    pm3.Normal('x', 0, 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(123)\n",
    "with model3:\n",
    "    hmc3 = pm3.HamiltonianMC(adapt_step_size=True)\n",
    "    \n",
    "point = {'x': np.array(.05)}\n",
    "trace3 = []"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "0.25\n",
      "10\n",
      "2.4338363922634105\n",
      "20\n",
      "0.8835295148661158\n",
      "30\n",
      "2.132369805866789\n",
      "40\n",
      "0.6638539397446028\n",
      "50\n",
      "0.7700032438774976\n",
      "60\n",
      "1.8083500351406778\n",
      "70\n",
      "2.01549699366806\n",
      "80\n",
      "0.620051066822726\n",
      "90\n",
      "1.2394151791390735\n",
      "100\n",
      "1.326608190951551\n",
      "110\n",
      "1.3458859508747132\n",
      "120\n",
      "1.4535020800208394\n",
      "130\n",
      "0.8470025039761001\n",
      "140\n",
      "0.9616054472149861\n",
      "150\n",
      "2.4700446424546842\n",
      "160\n",
      "1.579485048302955\n",
      "170\n",
      "2.126076524692899\n",
      "180\n",
      "1.4263912419505649\n",
      "190\n",
      "1.778245696911592\n",
      "200\n",
      "1.4859097351789183\n",
      "210\n",
      "0.9546259420166145\n",
      "220\n",
      "1.0072250233799054\n",
      "230\n",
      "0.9409202999644469\n",
      "240\n",
      "0.9783484856905771\n",
      "250\n",
      "1.0058739854796497\n",
      "260\n",
      "0.9099078292987351\n",
      "270\n",
      "1.5554539165456562\n",
      "280\n",
      "3.088081197599752\n",
      "290\n",
      "1.4375634311726395\n",
      "300\n",
      "1.1117815200764143\n",
      "310\n",
      "1.4225316718871723\n",
      "320\n",
      "2.2024230776071763\n",
      "330\n",
      "1.5205261622590753\n",
      "340\n",
      "1.3428281946994207\n",
      "350\n",
      "1.9597477496000024\n",
      "360\n",
      "1.665457343355618\n",
      "370\n",
      "0.9921134467269729\n",
      "380\n",
      "1.034038006686061\n",
      "390\n",
      "2.5815148520987616\n",
      "400\n",
      "0.8542902108760975\n",
      "410\n",
      "1.621225265014004\n",
      "420\n",
      "1.749029446162765\n",
      "430\n",
      "1.0936219616260738\n",
      "440\n",
      "1.103351819226634\n",
      "450\n",
      "1.2922387373823943\n",
      "460\n",
      "1.2931970508606732\n",
      "470\n",
      "1.6945914178093304\n",
      "480\n",
      "2.360210326424694\n",
      "490\n",
      "1.3767503315591088\n",
      "CPU times: user 152 ms, sys: 219 µs, total: 152 ms\n",
      "Wall time: 147 ms\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "for i in range(n_samples):\n",
    "    point, _ = hmc3.step(point)\n",
    "    trace3.append(point['x'])\n",
    "    if i % 10 == 0:\n",
    "        print(i)\n",
    "        print(hmc3.step_size)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x7fa496578b70>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAD8CAYAAABw1c+bAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt81PWd7/HXZyaZ3EMCuQAJIeFu5E4AlSoVlYK1UrfaCmpb3R7qVrfbbXu29vR0u3vc7Wm32+2pZ72stdVetFZFFAtFtOIFBSTI/R4ChJAQEkhICLlN5nP+SPSkMZBJSPKdyXyejwcPM/P7fSfv8UHe/PKd3+/7E1XFGGNM5PC4DmCMMWZgWfEbY0yEseI3xpgIY8VvjDERxorfGGMijBW/McZEGCt+Y4yJMFb8xhgTYaz4jTEmwkS5DtCVtLQ0zc3NdR3DGGPCxtatW6tUNT2YfUOy+HNzcyksLHQdwxhjwoaIHAt2X5vqMcaYCGPFb4wxEcaK3xhjIowVvzHGRBgrfmOMiTBW/MYYE2Gs+I0xJsJY8RtjTISx4jfGmAgTklfuGhOJntlc0uuxy+bm9GESM9jZEb8xxkQYK35jjIkwQRW/iCwSkQMiUiQiD1xkv9ki0ioit/Z0rDHGmIHRbfGLiBd4GFgM5ANLRST/Avv9GHi1p2ONMcYMnGCO+OcARaparKrNwLPAki72+1tgBXCqF2ONMcYMkGCKPws43uFxaftzHxGRLOAW4LGejjXGGDOwgil+6eI57fT4/wDfUdXWXoxt21FkuYgUikhhZWVlELGMMcb0RjDn8ZcCozo8zgbKOu1TADwrIgBpwI0i4g9yLACq+jjwOEBBQUGX/zgYY4y5dMEU/xZgvIjkASeA24FlHXdQ1bwPvxaRp4A/qupLIhLV3VhjjDEDq9viV1W/iNxP29k6XuBXqrpHRO5t3955Xr/bsX0T3RhjTG8EtWSDqq4B1nR6rsvCV9UvdzfWGPNxY0ue7/3gud/quyBm0LMrd40xJsJY8RtjTISx4jfGmAhjxW+MMRHGit8YYyKMFb8xxkQYK35jjIkwVvzGGBNh7J67xvSlwiddJzCmW3bEb4wxEcaO+I1xqLJR2H82iiNnmjlSFU2z+PB4oshJaCEvvpGkqIDriGYQsuI3poNnNpdc0vhl3otvV4UdJxs5crQYT3UxY/xFTJNTXC3nP9qnSaPYWTWG9wOTOBA3i5kjfExKbEC6uruFMb1gxW/MAKiqa2Dnnt0MOb2NaYH9TJcAjcRQlZBHfVIB3iGpVNS14NUWaK4np76Umc1/RJpfYfWRufwudjFL8iA5uvO9jozpOSt+Y/qLKsdKjnDm4EbyG7ezQFook+EcSl9IztjLSBiWTbbn//+KUHnkzF8ML/XXk161mYVntrCweSs/2n8XObnjyU9qGOh3YgYZK35j+pj6W9i/bwfxx99idOA4QzSBbYlXkz1xFtkjRjIyyDkbf1QC5cMXUJk2l6xjL/GDxl/xh+JPsi1nITNSG/v5XZjBzIrfmD7iazrNkb1rGXryXS6jjiKyeTP9LqZPncoVcdG9fl1/VALHxiylqeJNvnD6Tf54vJE93pu5PNnK3/ROUMUvIouAn9N2F60nVPVHnbYvAR4EAoAf+IaqbmjfdhSoA1oBv6oW9Fl6Y0KAr7GK4bsfZ1bli8RoM+/KDPy58/nEZTmM8/bRJ7Li4eTwBTRLLDdVvc7TRxM5MvaT5CU0983rm4jSbfGLiBd4GLiBtpunbxGRVaq6t8NufwZWqaqKyFTgOWBSh+3XqmpVH+Y2xrnYpiqy9j7OtJMriNYW1nquxjN+AdePTyG6n66QOZNxJcX+Bu6oeZ2Hj6Qw/LLJxHntlE/TM8Ec8c8BilS1GEBEngWWAB8Vv6qe67B/AqB9GdKYUOJpbWL0wSeZfvSX+LSJNXI1u8cuJ3v8VO6KfqN/v7kIlSMXQHMdX61/kR8ezeVTY+P793uaQSeY4s8Cjnd4XArM7byTiNwC/G8gA/h0h00KrBMRBf5LVR/vfVxjHFIlvfwNpu7+MZmt5awLzOaNUfcxKX8Go6MG8CJ4Ec7kLGLIwVLuPv9rVlX97cd/II25iGCKv6tJyo8d0avqSmCliFxD23z/9e2b5qlqmYhkAK+JyH5Vfftj30RkObAcICcnJ9j8xgyI5HPFTNr+r4yre5+DgSweSf1Xhs9YzNROH9pu7nRKZn9p9cZSlrOE/KO/Jqd8LdX1f0Nqgm9AvrcJf8EcppQCozo8zgbKLrRze6mPFZG09sdl7f89Baykbeqoq3GPq2qBqhakp6cHGd+Y/iXayrhDv2ThO7eSUbubn3ru5qW5zzHxqpsZcgln6vSF8wmj2J8yn8943uOlF37jNIsJL8EU/xZgvIjkiYgPuB1Y1XEHERkn0nZysojMBHzAaRFJEJGk9ucTgIXA7r58A8b0l4Tzx5n39p3MKfo/rG+dxvdHPUX69d8gOy3ZdbSP1I+4kgpPOvMP/zt7SipdxzFhotupHlX1i8j9wKu0nc75K1XdIyL3tm9/DPgc8EURaQEagC+0n+GTSdv0z4ff6xlVXdtP78WYPjG25HmSavaRe+IV/OrhQflvjB87mmUJu+FEaB23qCeK8hELmX7iaX773A/J/9Z/ILaoj+mGqIbeCTgFBQVaWFjoOoaJQL/feJjZG+9jXO1GtgfG8sqQO1kwSokK8S7NrnyblIqNbPr0Oq6bM911HOOAiGwN9jopW4/fmA+dP8Ocd/8b42o38rvAQrbm3M3CnNAvfYDhn/8Z0dJK42v/Smsg9A7mTGix4jcGCFQdpuahq8mu28G/yd0Mm/gJJg8Jn6tivWljKBvzeRY2/5n1m953HceEOCt+E/Gaj71P/aMLaG04yw+G/pgrJ+Uw1Od3HavHcpZ8n4B48b/xIzvqNxdlxW8iWv3e1wg8dRPVLdH8ac6vmTL3Bnye8CxNz5CRlI1fxg0t6/nzhg2u45gQZsVvItbZnauJfm4pRwKZ7Fn8And++rqwPyNm9M3fo1liiHr73wjYUb+5ACt+E5HOfPAS8S9+kUOaTfWtK1h85eA4E8aTlEHpuGXMb3mHTR/YmXGma1b8JuKc3r6GpFV/zT7NpXHZS1w1ZYLrSH0q96ZvExAvtesfch3FhCgrfhNRzuxdT8JLX6JIs5G7XmTWxFzXkfpcdEoWxSM+zfxza9lTVOw6jglBVvwmYpw5XEjMc0s5oWk0L1vBlHGjXUfqN9mf/gfipJmjf/q56ygmBFnxm4hwrqIY/d1t1Go8tbe9wLSJ41xH6lcJ2ZM5lDKPK6teoKzytOs4JsTYPXfN4FT45EdfNjXUc3b9wyQFGjl8+TeZ0/g+FHZ9kdPYkoFZVnkgpFz/bYa+cAtr1/yCkV96wHUcE0LsiN8MatrayrG3f0t6awXbxnyVOWPSXEcaMOmXX8tx3xjGHnmaxubwuyDN9B8rfjOo7Xj3FSa07OftzLuYf/ngndPvkgiNM77CeEp4/81XXKcxIcSmekxIemZzSa/HLpvbdge37dveZ3rtm7wZv5DrZg+O8/R7atx1X6Z28w/xFj4BC29xHceECDviN4NS0bES8kv/wAfeqVw1fxFhfkFur4kvgWOjP8fcpvfYt3+v6zgmRFjxm0GnurKM5F1PUinDyLt6KT5vZP81z138d3hQyv/8qOsoJkQENdUjIouAn9N2B64nVPVHnbYvoe0G6wHAD3xDVTcEM9aYrowteb5X4wKBACfWvch4PUfJ9G+SlRjXx8lC08WnxqKZEDOby0+9wlPvfB2fL+Yvtn44NWYiR7eHQiLiBR4GFgP5wFIRye+025+Baao6HbgHeKIHY43pM2ePbGVy6152Zt3O+FEjXMcJGUdH30qmVNN6cJ3rKCYEBPM78BygSFWLVbUZeBZY0nEHVT2n//8ejgmABjvWmL5SWXmKGxrWsjX2CmbPnOU6TkhpzruOKlKYVLbSdRQTAoIp/izgeIfHpe3P/QURuUVE9gOraTvqD3qsMZfqfGMzsyqe46SkcdknPus6Tujx+ihMWcwV/kLqK3t/xpQZHIIp/q7Oh/jYQt+qulJVJwGfpW2+P+ixACKyXEQKRaSwsrIyiFjGtFEF35HXSeMse0d+jvjYWNeRQlL1xC/gFWXIwRdcRzGOBVP8pcCoDo+zgbIL7ayqbwNjRSStJ2NV9XFVLVDVgvT09CBiGdOmrLSYeYFC3k74FMmpGa7jhCwdOpYdUVO54uwa/K12JW8kC+asni3AeBHJA04AtwPLOu4gIuOAw6qqIjIT8AGngZruxhpzKWrqG7n27Esc8OQRP3omAJuPDJ71doIV7FlQlcn5TDuzk8Stj5I1fHjbk3O/1Y/JTCjqtvhV1S8i9wOv0nZK5q9UdY+I3Nu+/THgc8AXRaQFaAC+0P5hb5dj++m9mAgTCCgpx9YSRxPlOTcR54ns8/WDkZAxhrOnE0ir3g7DF7mOYxwJ6jx+VV0DrOn03GMdvv4x8ONgxxrTF8pKi/mc7uTPSTeTmJjqOk5Y8Hij2BEzi7lN7/FOwwKS4nyuIxkH7BDJhKWzDc18snYVBySXhFFTXccJK4HMKcSIn/MnD7iOYhyx4jdhRxWijr5FMucpz/40YlM8PRKXPIwDkkd+/RYCgS5PsjODnP3EmLBTVnGSBYGNvBt/LXHJw1zHCUvHkmcxTkqpPGN354pEVvwmrDT5A8yoeoUTpBM7erbrOGErMXMcDeoj8fQu11GMA1b8JqxUH9tNnpSzP30x4o12HSdseaN9bI+axsyWD2hsOO86jhlgVvwmbFTXN3Jdw1q2eyaTmJHrOk7YO5s6mRSpZ9ebvVsJ1YQvK34TNmJK3saHn7OjrnMdZVAYkp5FlQ7Bs+s511HMALPiN2GhrOoM1wY28V7cfGITh7iOMyh4PB52xcxgSv1GzlSedB3HDCArfhPyAgElt+J1zmgSMTn2gW5f8qfn45NWDrzxG9dRzACy4jchr7z8BDPYzwfJ1+OJtitN+1LykGEc8+SQcuhF11HMALLiNyGtuVWZXr2WEjJJzLabt/U5EU7mLuEy/z5Kina7TmMGiBW/CWlVpYcYI2UcGHYD4vG6jjMojbnubgIqlL71lOsoZoBY8ZuQ1eQPMKfudQ5IHsmZea7jDFrpWWPZFzuVnNJX0EDAdRwzAKz4TciqLt1PllRxPH0+SFc3czN9pWHSbWTrSfYXvuE6ihkAQS3LbMxAa2xpZe65P7PXM5aktFHdDzC9V/gkk9KiadRozr75f8FzvPsxHyq4u/9ymX5jR/wmJJ0t3ctwqaYs45N2tD8AEuNi2BM7g0n1hTT77baMg11QxS8ii0TkgIgUicgDXWy/Q0R2tv95T0Smddh2VER2ich2ESnsy/BmcDpXX88V9evZLRNISstyHSdiRI0qIEXOse9gkesopp91W/wi4gUeBhYD+cBSEel8Xt0RYL6qTgUeBB7vtP1aVZ2uqgV9kNkMcttWPUKG1HAy4xOuo0SUy8aNoUYTaTmx3XUU08+CmeOfAxSpajGAiDwLLAH2friDqr7XYf9NQHZfhjSRo6m5ibwDv+CA5JI0zI72B0LHm9M3Rk+joHEr7x5aQFRU9/Uw1w7lwlIwUz1ZQMdPe0rbn7uQvwb+1OGxAutEZKuILO95RBNJtq55kmwqODb0apvbd+Bc6iQSpJGaUz34gNeEnWCO+Lv66evyfm0ici1txd/xd/R5qlomIhnAayKyX1Xf7mLscmA5QE5OThCxzGDT2hogc8cjlHhzSLZll51ITRvB6VPJDKvdAyPt2onBKpgj/lKg4/l02UBZ551EZCrwBLBEVT+6n5uqlrX/9xSwkrapo49R1cdVtUBVC9LT04N/B2bQ2PrnFxirx6ie8TXEY0f7Lng9Hnb6pjPNv5vmlhbXcUw/Cab4twDjRSRPRHzA7cCqjjuISA7wInCXqh7s8HyCiCR9+DWwELAFQUyXfFseoUpSmfype1xHiWgNqZOIk2ZqTx1zHcX0k26LX1X9wP3Aq8A+4DlV3SMi94rIve27/SMwDHik02mbmcAGEdkBvA+sVtW1ff4uTNjbt2MT01u2cXTsnXijY1zHiWipwzKp0FTSa/e4jmL6SVBX7qrqGmBNp+ce6/D1V4CvdDGuGJjW+XljOqt54+c0qI9JN/2d6ygRz+MRdsdM5+qmd9jcvAifz/4hHmzsyl3jXEVZCTNrXmNPxk0kptjnO6GgaegkfOKn9tRR11FMP7DiN84dXPMQMdJC1qK/dx3FtEtNTeOEpjHcpnsGJSt+41RjYyMTS19gT9xsRoyd6jqOaefxCHtjpzM1sI/GpkbXcUwfs+I3Tm1//fdkUI1nzsc+IjKO+YdOIlpaqa844jqK6WNW/MapuB1PUSHpTLrmVtdRTCepKamUaCYjz9kZ2IONFb9xpmjvNqa1bOf4mC8gXrs1RKgRj7AvdjpTAvtpbDzvOo7pQ1b8xpmTbzxCs3qZsOhrrqOYC2gdNgGvKPV2MdegYodZpv8UPnnBTXWNLUyp/CN74mYx49g6sF4JSUOHpFB6Ip3h5/bSwmWu45g+Ykf8xold+/YxRM6TOGau6yjmIsQj7IuZxuTAfpqbm1zHMX3Eit84kVi+kXLSGZdrK0CGupah4/FJK3WnSlxHMX3EpnpMv+l4g4+OTtfWcWNgH6/H30jJsZoBTmV6KjU1jZPlQ0mr24sy3nUc0wfsiN8MOE/FHlpV8A23OeNw4PEIe31TmNq615ZqHiSs+M2AamlVZja9zw7vZGLiElzHMUFqSJ1IrLRQW2l35hoMrPjNgKo6dYIMqaEydbrrKKYHUodlcFqTGVq733UU0wes+M2ASq/ZQY0mMiTDbq8ZTrweD7uipzLVv5sWf6vrOOYSWfGbAXOusZmC1p3siJmFeLyu45geqk+ZQII0crbyhOso5hIFVfwiskhEDohIkYg80MX2O0RkZ/uf90RkWrBjTeQ4V3GYGGnBn57vOorphdS0EZzVeJJtuifsdVv8IuIFHgYWA/nAUhHp/JN7BJivqlOBB4HHezDWRIgx57ZxlBEkDklzHcX0gtfrZVfUFKa27KI1YNM94SyYI/45QJGqFqtqM/AssKTjDqr6nqpWtz/cBGQHO9ZEhtO155hCEYcSZoGI6ziml84OmcQQqae6qtx1FHMJgin+LKDjOVyl7c9dyF8Df+rlWDNIyam9BFSIyZjgOoq5BEPSs6nXGBKqD7qOYi5BMMXf1eGZdrmjyLW0Ff93ejF2uYgUikhhZWVlELFMuAgElPzGbez2TCAmPtF1HHMJoqO87IyazJSWnQQCAddxTC8FU/ylwKgOj7OBss47ichU4Algiaqe7slYAFV9XFULVLUgPd1uuD2YVFZXM1oqqEie4jqK6QNnki5jmNRSc/qU6yiml4Ip/i3AeBHJExEfcDuwquMOIpIDvAjcpaoHezLWDH6+0/tpUS/xGWNcRzF9IDljFI0aTWyNTfeEq24XaVNVv4jcD7wKeIFfqeoeEbm3fftjwD8Cw4BHpO2DO3/70XuXY/vpvZgQ5A8oU5u3scubT5Qv1nUc0wd80dHs9OZzWdNONNBq12SEoaBW51TVNcCaTs891uHrrwBd3i27q7EmclRWVTFPTrN7yA0McR3G9JlTifnMqd1B0fYNjJs533Uc00N25a7pVwnV+2jUaBLTc11HMX0oMSOHFvVyuvB511FML1jxm37T0qrMaNnGrqjL8Ub7XMcxfSguJoadnklkl78G2uWJeiaEWfGbflNVVUGanKVmyOWuo5h+UJaQT5aepGTf+66jmB6y4jf9Jrl6L+c1hkRbiXNQik/Po1WFk5uecx3F9JAVv+kXDQ0NzPBvZ2f0FDzeaNdxTD9IjI9ln28ymSfWuY5iesiK3/SLXe+8TKqcozbF1uQbzM7mLmZ0awnlh3e6jmJ6wIrf9IvAzheo1XiS07K739mErdGf+AIApe/9wXES0xNW/KbP1Z2rY3LdBnb5piHeoC4VMWEqe/Q49nknMuzYWtdRTA9Y8Zs+t+etFSRKAw2pk1xHMQOgatQixviLqCo94DqKCZIVv+lznj0vUk0yicNsBe5IkHXV5wE4tsGme8KFFb/pUzU1Z5hSv5HijBsQj/31igR54y/nkCeP5OI/db+zCQn2k2n61N43nyNOmkmdc7vrKGaAiAjlIxcyvnkvtRUlruOYIFjxmz4Vs38llTKMvJnXuY5iBlD63NsAKN7wrOMkJhhW/KbPVFaeZErDFkqGL7SleiPMxMsLOCLZxB5a7TqKCYIVv+kzB998Fp+0kn7lMtdRzADzeISSzOsZ37CD89UnXccx3bDiN30m4dDLlHsyyZlytesoxoGUWbfiFeXwO7Z2T6gLqvhFZJGIHBCRIhF5oIvtk0Rko4g0ici3O207KiK7RGS7iBT2VXATWk6WlTC5aTsnsm6EtruwmQhz+YyrOE4m3gOvuI5iutFt8YuIF3gYWAzkA0tFpPMCLGeArwP/foGXuVZVp6tqwaWENaGr6M1niJIAI+bd4TqKcSQqyktx2gLGn9tKY22V6zjmIoI54p8DFKlqsao2A88CSzruoKqnVHUL0NIPGU0YSCl+hRLvKLIm2r/tkWzI7NuJllYOvfm06yjmIoIp/izgeIfHpe3PBUuBdSKyVUSW9yScCQ/HjxaR37KHUzmftmmeCDel4BqOMRLfvhddRzEXEUzxd/WT3JN7rc1T1Zm0TRXdJyLXdPlNRJaLSKGIFFZWVvbg5Y1rR995Bo8oo66+03UU45jX66F4+CLGn99BfdXx7gcYJ4Ip/lJgVIfH2UBZsN9AVcva/3sKWEnb1FFX+z2uqgWqWpCenh7sy5sQkHb0jxRHjSVzzBTXUUwISL9iGR5RDtt0T8gKpvi3AONFJE9EfMDtwKpgXlxEEkQk6cOvgYXA7t6GNaGn+NAeLms9QHXeTa6jmBCRP7WAg5JHwsGVrqOYC+i2+FXVD9wPvArsA55T1T0icq+I3AsgIsNFpBT4JvA/RaRURJKBTGCDiOwA3gdWq6ot3D2IlL3zGwBGz7dpHtPG4xFKs29kbPN+assOuo5juhDUXTJUdQ2wptNzj3X4+iRtU0Cd1QLTLiWgCV0aCDDq+B/ZHzOZSdkTXMcxIWTkVXfAHx7m6Fu/ZerSB13HMZ3Ylbum1/Zte5fRWsr5CZ9zHcWEmImT8tnlmURKsV3MFYqs+E2vVW/6Hc3qZcJ1d7mOYkKMiHBq9E3ktBzhzJEdruOYTqz4Ta+0tLQwofJV9iddQWKKnYVlPm701XfQqkLpO791HcV0YsVvemXXhldIpxqZ+nnXUUyIGjdmDDuip5Jx7BXQnlz6Y/qbFb/pleZtf6COOCZdc5vrKCaEncm7meGtJ6k8sNF1FNNBUGf1GNPRuXN1TD77FgeGLWBWbILrOMalwicvunl8RiJNB6MoX/vvpJ/7zF9uLLi7H4OZi7EjftNju9f/gURpIGn2UtdRTIgbnRrDB95pjKrZgra2uo5j2lnxmx7z7XmeShnK+DmLXUcxYaBx5FxSqaXkyAHXUUw7m+oxPVJZUcaUhi3syLqddK/99Yl0m4+c6XYfjR9GlSZTXbSJk97hHz0/11bwdsaO+E2PHHjjt0RLKxnzvug6igkTCdHCpug5XN6yB2mudx3HYMVvemjo4Zco8eaQkz/XdRQTRuqHTSFaWmmssOmeUGDFb4J27NBO8v17qchdYjdcMT2SMyyJ3ZrH6Lptdk5/CLDiN0E78cYTtKqQd/1XXEcxYcYjsDduDqP1BK11Fa7jRDz7dM5cXPt52s3+ViaUvcRu3zSmnVgPJxznMmEnYeQ46g/HEF2xnUDyItdxIpoVv7moD8/aOFV+jM/IWbYkT6MxiDM5jOksLc7L2965XNu0kZ3++a7jRDSb6jFByazeRqWmkJIxqvudjbmAmmEziZUWmk/ucx0logVV/CKySEQOiEiRiDzQxfZJIrJRRJpE5Ns9GWtCX93588wK7GJH3Gw8HjtWML2Xk57MDh1H3tkt9iGvQ93+FIuIF3gYWAzkA0tFJL/TbmeArwP/3ouxJsT5y/fiFUUyLncdxYS5KIE9CVeQRQWndq5zHSdiBXP4NgcoUtViVW0GngWWdNxBVU+p6hagpadjTWhrDbQys+E9PpDLSUxKdh3HDAKpI/I4rUnUrH/IdZSIFUzxZwHHOzwubX8uGJcy1oSA6vKjZEgN5UNnu45iBomhscKb0dcwoWYD9WU21+9CMMXf1ZU6wU7OBT1WRJaLSKGIFFZWVgb58qa/ja7ZzHHNYKh9qGv60oipNGk0Jat/6jpJRAqm+EuBjj/12UBZkK8f9FhVfVxVC1S1ID3dbuUXCop2bGAyRexOvBKPx67UNX0nOzmaDfELGHPiZVrPVbmOE3GCKf4twHgRyRMRH3A7sCrI17+Uscax6vX/yXmNIXHkJNdRzCDku/rrxNDM4TU/dx0l4nRb/KrqB+4HXgX2Ac+p6h4RuVdE7gUQkeEiUgp8E/ifIlIqIskXGttfb8b0ndNlR5hWvY7Nvrn4fDGu45hB6Kor5rHJM5OMfb9GbdXOARXUlbuqugZY0+m5xzp8fZK2aZygxprQd/iVf2MmSvOIAuJdhzGDktcj1M7+Oimbv8zhtf/J2Ju/4zpSxLCrcczH1J6pYHLZCj5IXkBKUqLrOGYQm3/DzRTKZIZufxRtPu86TsSw4jcfs//lnxIvTaQutCMw079iorycLvh7UgPVHFn3iOs4EcOK3/yFhnNnmXDsabbFXcH4KXNcxzERYP7CW/hA8kn54BFoaXQdJyJY8Zu/sHflj0nhHL5r/8F1FBMhYqO9VMz4O4YGTnNkrV3NOxCs+M1HzlWfYuLhJymMvZLL51znOo6JINcuvo33PdMY9sFDBM5Xu44z6Fnxm4/sf/6fiNcGEm/8X66jmAgTG+2l7urvkxg4x+GV/+I6zqBnxW8AOH2imCknnmNz8kImTbW5fTPwrp1/PW/EfJKcQ7+h6fQx13EGNSt+A8DxFd8FlKzP2tG+ccPjEZIW/xOoUvL8d13HGdSs+A0l215n+pm1bMxcSs5YW57BuDN3xnQY3FGLAAAOYUlEQVReHXIr40+upnrvm67jDFpW/BFO/c2w+luUk8a0pQ+6jmMM05Y+yAlNo+Hlv4dWv+s4g5IVf4Tb+eJPyPEf5fCs75Oamuo6jjGMHpFO4WXfYWRTMUfW/IfrOIOSFX8Eqy07xLi9D7HVN5urbvyi6zjGfORTf3UPGz0zydz6Uxrtg94+F9QibSa8bX7+4ze70ECAIQeeJVuF8sxP4n/xZw6SGdO1WF8Unpt+ir68iLLfLGfMN9aC2D0h+ood8UeoptJtXBYo4k+Jt5CeHOc6jjEfM3fmTF7Luo8xZzdx7PXHuh9ggmbFH4Ga6qq4svZVNjCDUTljXccx5oKuu+sBtspk0t79ZxqrbMqnr9hUT4TRlgbGlTxPNUmcy/sUqXZLReNIV1OQXTky/EYmlj1EyaO3UD1xKYiHubd9q5/TDW5BHfGLyCIROSAiRSLyQBfbRUQeat++U0Rmdth2VER2ich2ESnsy/Cmh1SJP7yGND3D+rRlpMb7XCcyplt5Q2NZEf95JrQehpLNruMMCt0Wv4h4gYeBxUA+sFRE8jvtthgY3/5nOfBop+3Xqup0VS249Mimt1qObWJq6x5Wxn2OvOHDXMcxJmjjcnN5lasoqHud1upS13HCXjBH/HOAIlUtVtVm4FlgSad9lgC/0TabgBQRGdHHWc0laCrfyyfqX+N1mUdO3kTXcYzpkSgPtI5ZwDGGM7FsBY1nrPwvRTDFnwUc7/C4tP25YPdRYJ2IbBWR5b0NanqvseooV51+kc1cjmf8Ajw2r2/C0LA4D5szlxKrjZz8xW2o3bSl14Ip/q5aQnuwzzxVnUnbdNB9InJNl99EZLmIFIpIYWVlZRCxTDBO7lrP7JN/oIhR1Iz9LAnRVvomfI1JT2RF0h3kNuzlwK++Ctq5ikwwgin+UmBUh8fZQFmw+6jqh/89BaykberoY1T1cVUtUNWC9PT04NKbiyrf9SZJK26nglQOj76dlFiv60jGXLIJOdmsGXonk8pf4sBLP3IdJywFU/xbgPEikiciPuB2YFWnfVYBX2w/u+cK4KyqlotIgogkAYhIArAQ2N2H+c0FnNj1JskrbqeKVHbmfIm0pBjXkYzpEyLwya/+jHd985i440cUv/U715HCTrfFr6p+4H7gVWAf8Jyq7hGRe0Xk3vbd1gDFQBHwC+Br7c9nAhtEZAfwPrBaVdf28XswnRzZ+BJDV3ye0wzBf+cqhidb6ZvBJT7Gx8S/eYZdnklkrf8GJ3a84TpSWAnqAi5VXUNbuXd87rEOXytwXxfjioFpl5jR9MDeV3/B+Pe+Q7FnNL4vrWBs7hiqtrlOZUzfS0tN4fzdL1D+y4UMW3kHFbEryZx4hetYYcGWbBgkNBBgx+8eIH/jt9kbnc/Qr60jL3eM61jG9KucUaNoXraCOuKJ+f1tVBbvcB0pLFjxDwLNDfXs/PmtTCt6lA0JNzDmG2uxD8hNpJgwIZ+a21bQoh7kt0uoPLLLdaSQZ2v1hLnyI3tofPpOpvmLeSvnPj7xpX/B67V/z83g1tU6PwdHfpl5ZU/ifepG3sq5i9jkrq9Ot3V+7Ig/rG1b91sSfn09Q1sq2HLVfzH/nh9a6ZuINWJoMhuz7qYFL9NKfk39Wbse6EKsJcJQXW0Nm39+FzPeu5+KqJGc+/J6Zi+83XUsY5zLTE2mcNTdnCeGucefpK7KlnboihV/mNm5cR1nfnYVs8+8wuYRd5Hz3zeQZWvvGPORtCGJ7Mq9hzLSuObkr6ktO+Q6Usix4g8Tp8+c4a2HvsLktZ8nhmaKFv2OuV/9T2Ji7O5ZxnSWkhjH8XF3slfGcUP17zlfvAkN2PIOH7LiD3Et/lbeWPFftDxUwPwzz7NjxK2kfKuQCVfe5DqaMSEtIdZH3cRbedN7Fdc2rMNzcDV+f4vrWCFBNAQXOSooKNDCwsi+Z0sgoLz79jqS3vlfTG/dTTHZHBu5mPhUW+3amJ7QgFJ5dCc3nl/FcRlO9LLfkzVhhutYfU5EtgZ7zxM7nTPEtAaUtza8RdTbP+Ya/3vUSDL7Z/6AmoZW4j32C5oxPSUeIWPMNFZXZDCv8g/EPf0pts36J2Z85m/aFv6JQNYkIaKpxc/611fz3g8XseCNJRT4t3Fg0tdI+ofdTLr5m4iVvjGXJDNzBNvGLKfYN4EZH3yXbT+9mbNV5a5jOWFTPQOoq4tOKhsCnCsvYtr597hMSqgjnj2JV+EdOZ1AdLyDlMYMbrM++3U2P/PPzC5+hDpJ4tjs7zNj8T1hf3BlUz0hrtEPZVU1JFfvYn7rJpKkgaOebN4ecgsxmROQqBgCrkMaM0hFRUcz70v/QtHOTxN4+W+ZueVb7Nr1DCl/9R+MmjDddbwBYcU/QBqa/ByuOEtCzX5mtnzAfKmkWaPYFTOd5ozpeJJHEBOh843GuDBu6pX4L9vIxud/wpQDDxH79LW8n76EMbc+SNrwUd2/QBiz4u9HZVVnOLh5La0HX2Nczbsskwpa1Mv+qIkcHTKf2PQ8AlFx9kGLMY5ERUdz5bL/QVXFF9n9/PeZVfky/kdXsznjZkZ/5jsMz5ngOmK/sOLvQ2frGzi4411q9r1FSvk7TGnZzSelhSZ8HEuewTvR84hNG0MgOh4f2HSOMSEiLTObtPufpOTQTk6u/t/MPLUS+eWLbEu8Eu/se5j8ic/iiRo8dRnUh7sisgj4OeAFnlDVH3XaLu3bbwTOA19W1Q+CGduVcPhwNxBQSkqPU35gMw2HNzKkaisTW/aTKI0AlEVlc3rEfFKnLiZr2nWIL77LD3eNMQMrmNU5Tx4v4sianzGxfBVDqaWKFA6nLSBh+i2ML7iemNjQO/GiTz/cFREv8DBwA203Vd8iIqtUdW+H3RYD49v/zAUeBeYGOTak+f2tnDx5gtPHD1BXfgh/xX4Sq/eR3VRErpwhFwioUOrLo3jkZ4gZO4+c6dcxMi2Hka7DG2N6ZfiocQz/6sM0Nf6ErW88i+x5kWmVrxD7+os0vOZjV+xkzmXMIj63gMwJBWSMzMPj9bqOHbRgfneZAxS130YREXkWWAJ0LO8lwG/ab8G4SURSRGQEkBvE2D4VaA3QGvDT6vfj97fQ6vfT6m+m1d9CIOCntaUFv7+FpoZ6ms/X0nK+lpaGOvwNtbSeq0TqK/E2nia26QyJzVVkBirIlgay21/frx7KonOoGDabk5lTGJI3k1H5V5KTOJSc/npTxhgnYmLjmXXjPXDjPdTX1bB90xoaD75B5un3ubzkCTzHfwHvQLNGcdKbSXVMFo2JOWhCJhKfSlTiUHyJacQlpRIVG483OoZoXyxRvhiifXFE+3z4YuLweLyIyICdUhpM8WcBxzs8LqXtqL67fbKCHNtnnvrBMr4sq/EA0b18jSaiqZEU6rwpNMYP52DSHDxD84gbPo6hWRNIGzWRHF+clbwxESYhKYXpNyyDG5YBUF9Xw7HdG6k7vgc9U4yvroQhjaWMPbWLRGno1feoIoW0fzrWl7G7FEzxd3WOYecPBi60TzBj215AZDmwvP1hk4jsDiLbx9zdm0Efc7qnA9KAqj751gPHMve/cMsL4Ze5F3m/3S9BeuAimWvhn3t9WvfoYHcMpvhLgY4ntWYDZUHu4wtiLACq+jjwOICIFAb7IUUoCLe8YJkHQrjlhfDLHG55ITQyBzOhtAUYLyJ5IuIDbgdWddpnFfBFaXMFcFZVy4Mca4wxZgB1e8Svqn4RuR94lbZTMn+lqntE5N727Y8Ba2g7lbOIttM5777Y2H55J8YYY4IS1BUJqrqGtnLv+NxjHb5W4L5gxwbh8R7u71q45QXLPBDCLS+EX+ZwywshkDkkV+c0xhjTf2yZGGOMiTAhXfwi8m0RURFJc52lOyLyoIjsFJHtIrJOREL+wl0R+YmI7G/PvVJEUlxnuhgRuU1E9ohIQERC+kwOEVkkIgdEpEhEHnCdpzsi8isROdXb06gHmoiMEpH1IrKv/e/E37nO1B0RiRWR90VkR3vmf3aVJWSLX0RG0bbUQ4nrLEH6iapOVdXpwB+Bf3QdKAivAZNVdSpwEPiu4zzd2Q38FfC26yAX02GpksVAPrBURPLdpurWU8Ai1yF6wA98S1UvA64A7guD/8dNwAJVnQZMBxa1nwU54EK2+IGfAf/ABS74CjWqWtvhYQJhkFtV16mqv/3hJvhoZYqQpKr7VPWA6xxB+GiZE1VtBj5cqiRkqerbwBnXOYKlquUfLgSpqnXAPtpWCghZ2uZc+8Po9j9OeiIki19EbgZOqOoO11l6QkT+VUSOA3cQHkf8Hd0D/Ml1iEHiQkuYmH4gIrnADGCz2yTdExGviGwHTgGvqaqTzM4WmBaR14HhXWz6HvA/gIUDm6h7F8usqi+r6veA74nId4H7gR8MaMAudJe5fZ/v0far89MDma0rweQNA0EvVWIujYgkAiuAb3T6rTskqWorML3987SVIjJZVQf8cxVnxa+q13f1vIhMAfKAHW3L/JMNfCAic1T15ABG/JgLZe7CM8BqQqD4u8ssIl8CbgKu0xA4t7cH/49DWTDLnJhLJCLRtJX+06r6ous8PaGqNSLyJm2fqwx48YfcVI+q7lLVDFXNVdVc2n6IZrou/e6IyPgOD28G9rvKEqz2m+R8B7hZVc+7zjOI2FIl/az95k+/BPap6n+4zhMMEUn/8Mw5EYkDrsdRT4Rc8YexH4nIbhHZSds0VcifXgb8J5AEvNZ+Gupj3Q1wSURuEZFS4EpgtYi86jpTV9o/MP9wqZJ9wHOhvlSJiPwe2AhMFJFSEflr15m6MQ+4C1jQ/nd3u4jc6DpUN0YA69s7Ygttc/x/dBHErtw1xpgIY0f8xhgTYaz4jTEmwljxG2NMhLHiN8aYCGPFb4wxEcaK3xhjIowVvzHGRBgrfmOMiTD/D3xOn0+dXJxaAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import seaborn as sns\n",
    "sns.distplot(trace)\n",
    "sns.distplot(trace3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There still seems to be a problem here where in the PyMC4 implementation, the step_size keeps getting smaller and smaller, causing the sampler to take very long. Haven't figured it out yet."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.7264840418426366"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hmc.step_size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.7268868356485505"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hmc3.step_size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.97779363])"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hmc.potential._stds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.97756585])"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hmc3.potential._stds"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
